!MSc Thesis Code
!Author: Amin Deldari Alamdari
!Topic: Optimization of rib reinforced I-beams
!Department of Mechanical Engineering, Bogazici University

!-------------------------------------------------------------
!Units: m, kg, m/s^2, kg/m^3, N, Pa
!-------------------------------------------------------------!------------------------------------------
FINISH
/CLEAR,NOSTART

!*******************************************************************************************************
!								Initial OF THE PROGRAM
!*******************************************************************************************************
*CREATE,Initial

KEYW,PR_SGVOF,0 			!NO WARNING MESSAGE
/NERR,0,100000,,0,0 		!LIMITING THE NUMBER OF WARNING AND ERROR MESSAGES DISPLAYED
!/FTYPE,ALL,INT 
/PREP7
/NOPR   
OUTRES,ALL,NONE			!CONTROLLING THE SOLUTION DATA WRITTEN TO THE DATABASE
OUTRES,MISC
OUTRES,ESOL
OUTPR,ALL,LAST 
!/FILNAME,buckling,0
KEYW,PR_SET,1   
KEYW,PR_STRUC,1 
KEYW,PR_THERM,0 
KEYW,PR_FLUID,0 
KEYW,PR_ELMAG,0 
KEYW,MAGNOD,0   
KEYW,MAGEDG,0   
KEYW,MAGHFE,0   
KEYW,MAGELC,0   
KEYW,PR_MULTI,0 
!*
/VIEW,1,1,1,1   
/ANG,1  
/REP,FAST  
!*
*AFUN,DEG

!*******************************************************************************************************
!										ELEMENT TYPE
!*******************************************************************************************************


ET,1,SHELL181				!ELEMENT TYPE

KEYOPT,1,1,0 
KEYOPT,1,3,0 
KEYOPT,1,5,0 
KEYOPT,1,8,0 
KEYOPT,1,9,0 

!*******************************************************************************************************
!										MATERIAL PROPERTIES
!*******************************************************************************************************

ThicknessPly=0.000125				!PLY THICKNESS
LaminatePlyNumberTop=8				!PLY NUMBER OF THE TOP FLANGE LAMINATES
LaminatePlyNumberWeb=4				!PLY NUMBER OF THE WEB LAMINATES
LaminatePlyNumberBottom=8			!PLY NUMBER OF THE BOTTOM FLANGE LAMINATES

ElasticModulus1=165e9				!LONGITUDINAL ELASTIC MODULUS (E1) (FIBER DIRECTION)
ElasticModulus2=9e9					!TRANSVERSE ELASTIC MODULUS (E2)
ElasticModulus3=9e9					!THROUGH-THICKNESS ELASTIC MODULUS (E3)

ShearModulus12=5.6e9				!IN-PLANE SHEAR MODULUS (G12)
ShearModulus13=5.6e9				!TRANSVERSE SHEAR MODULUS (G13)
ShearModulus23=2.8e9				!THROUGH-THICKNESS SHEAR MODULUS (G23)

PoissonsRatio12=0.34				!MAJOR POISSON'S RATIO (v12)
PoissonsRatio13=0.34				!MAJOR TRANSVERSE POISSON'S RATIO (v13)
PoissonsRatio23=0.5					!THROUGH-THICKNESS POISSON'S RATIO (v23)

TensileStrength1=2560e6				!LONGITUDINAL TENSILE STRENGTH (XT) (FIBER DIRECTION)
CompressiveStrength1=-1590e6		!LONGITUDINAL COMPRESSIVE STRENGTH (XC)
TensileStrength2=73e6				!TRANSVERSE TENSILE STRENGTH (YT)
CompressiveStrength2=-185e6			!TRANSVERSE COMPRESSIVE STRENGTH (YC)
TensileStrength3=63e6				!THROUGH-THICKNESS TENSILE STRENGTH (ZT)
CompressiveStrength3=-185e6			!THROUGH-THICKNESS COMPRESSIVE STRENGTH (ZC)

ShearStrength12=90e6				!IN-PLANE SHEAR STRENGTH (S12)
ShearStrength13=90e6				!TRANSVERSE SHEAR STRENGTH (S13)
ShearStrength23=57e6				!THROUGH-THICKNESS SHEAR STRENGTH (S23)

TensileStrain1=0.01551				!LONGITUDINAL TENSILE FAILURE STRAIN (E1T) (FIBER DIRECTION)
CompressiveStrain1=-0.011			!LONGITUDINAL COMPRESSIVE FAILURE STRAIN (E1C)
TensileStrain2=0.0081				!TRANSVERSE TENSILE FAILURE STRAIN (E2T)
CompressiveStrain2=-0.032			!TRANSVERSE COMPRESSIVE FAILURE STRAIN (E2C)
TensileStrain3=0.007				!THROUGH-THICKNESS TENSILE FAILURE STRAIN (E3T)
CompressiveStrain3=-0.032			!THROUGH-THICKNESS COMPRESSIVE FAILURE STRAIN (E3C)

ShearStrain12=0.05					!IN-PLANE SHEAR FAILURE STRAIN (Gama12)
ShearStrain13=0.05					!TRANSVERSE SHEAR FAILURE STRAIN (Gama13)
ShearStrain23=0.021					!THROUGH-THICKNESS SHEAR FAILURE STRAIN (Gama23)

Density=1600						!MATERIAL DENSITY

Theta1=0
Theta2=90
Theta3=45
Theta4=-45

MPTEMP,,,,,,,,   
MPTEMP,1,0   
MPDATA,EX,1,,ElasticModulus1 
MPDATA,EY,1,,ElasticModulus2 
MPDATA,EZ,1,,ElasticModulus3 
MPDATA,PRXY,1,,PoissonsRatio12  
MPDATA,PRYZ,1,,PoissonsRatio23 
MPDATA,PRXZ,1,,PoissonsRatio13 
MPDATA,GXY,1,,ShearModulus12
MPDATA,GYZ,1,,ShearModulus23
MPDATA,GXZ,1,,ShearModulus13
MPDATA,DENS,1,,Density

TB,FCLI,1,1,9,1
TBTEMP,0
TBDATA,1,TensileStrength1
TBDATA,2,CompressiveStrength1
TBDATA,3,TensileStrength2
TBDATA,4,CompressiveStrength2
TBDATA,5,TensileStrength3
TBDATA,6,CompressiveStrength3
TBDATA,7,ShearStrength12 
TBDATA,8,ShearStrength23
TBDATA,9,ShearStrength13

TB,FCLI,1,1,9,2
TBTEMP,0
TBDATA,1,TensileStrain1
TBDATA,2,CompressiveStrain1
TBDATA,3,TensileStrain2
TBDATA,4,CompressiveStrain2
TBDATA,5,TensileStrain3
TBDATA,6,CompressiveStrain3
TBDATA,7,ShearStrain12 
TBDATA,8,ShearStrain23
TBDATA,9,ShearStrain13

FC,1,S,XTEN,TensileStrength1
FC,1,S,YTEN,TensileStrength2
FC,1,S,ZTEN,TensileStrength3
FC,1,S,XCMP,CompressiveStrength1
FC,1,S,YCMP,CompressiveStrength2
FC,1,S,ZCMP,CompressiveStrength3
FC,1,S,XY,ShearStrength12
FC,1,S,YZ,ShearStrength23
FC,1,S,XZ,ShearStrength13

FC,1,EPEL,XTEN,TensileStrain1
FC,1,EPEL,YTEN,TensileStrain2
FC,1,EPEL,ZTEN,TensileStrain3
FC,1,EPEL,XCMP,CompressiveStrain1
FC,1,EPEL,YCMP,CompressiveStrain2
FC,1,EPEL,ZCMP,CompressiveStrain3
FC,1,EPEL,XY,ShearStrain12
FC,1,EPEL,YZ,ShearStrain23
FC,1,EPEL,XZ,ShearStrain13


!*******************************************************************************************************
!								BEAM THICKNESS AND PLY NUMBER PROPERTIES
!*******************************************************************************************************

ThicknessTop=(LaminatePlyNumberTop*ThicknessPly)*12				!TOP THICKNESS (1 mm)*n
ThicknessWeb=(LaminatePlyNumberWeb*ThicknessPly)*4				!WEB THICKNESS (0,5 mm)*n
ThicknessBottom=(LaminatePlyNumberBottom*ThicknessPly)*12		!BOTTOM THICKNESS (1 mm)*n


ThicknessFlangeTotal=ThicknessTop+ThicknessBottom				!TOTAL FLANGE THICKNESS

PlyNumberTop=ThicknessTop/ThicknessPly							!TOTAL PLY NUMBER OF THE TOP FLANGE
PlyNumberWeb=ThicknessWeb/ThicknessPly							!TOTAL PLY NUMBER OF THE WEB
PlyNumberBottom=ThicknessBottom/ThicknessPly					!TOTAL PLY NUMBER OF THE BOTTOM FLANGE

!*******************************************************************************************************
!										SHELL THICKNESSES
!*******************************************************************************************************

SECT,1,shell,,Top Flange 
*DO,ply,1,PlyNumberTop/LaminatePlyNumberTop,1
SECDATA, ThicknessPly,1,Theta3,3    
SECDATA, ThicknessPly,1,Theta1,3   
SECDATA, ThicknessPly,1,Theta4,3 
SECDATA, ThicknessPly,1,Theta1,3    
SECDATA, ThicknessPly,1,Theta1,3    
SECDATA, ThicknessPly,1,Theta4,3 
SECDATA, ThicknessPly,1,Theta1,3   
SECDATA, ThicknessPly,1,Theta3,3    
SECOFFSET,TOP    
SECCONTROL,,,, , , , 
*ENDDO

SECT,2,shell,,Web
*DO,ply,1,PlyNumberWeb/LaminatePlyNumberWeb,1
SECDATA, ThicknessPly,1,Theta3,3    
SECDATA, ThicknessPly,1,Theta4,3   
SECDATA, ThicknessPly,1,Theta4,3   
SECDATA, ThicknessPly,1,Theta3,3    
SECOFFSET,MID    
SECCONTROL,,,, , , ,
*ENDDO

SECT,3,shell,,Bottom Flange  
*DO,ply,1,PlyNumberBottom/LaminatePlyNumberBottom,1
SECDATA, ThicknessPly,1,Theta3,3    
SECDATA, ThicknessPly,1,Theta1,3   
SECDATA, ThicknessPly,1,Theta4,3 
SECDATA, ThicknessPly,1,Theta1,3 
SECDATA, ThicknessPly,1,Theta1,3    
SECDATA, ThicknessPly,1,Theta4,3 
SECDATA, ThicknessPly,1,Theta1,3   
SECDATA, ThicknessPly,1,Theta3,3    
SECOFFSET,BOT    
SECCONTROL,,,, , , , 
*ENDDO

*END

!*******************************************************************************************************
!										STARTING PARAMETERS
!*******************************************************************************************************

*CREATE,StartingPar

n=4        					!DIMENSION OF THE PROBLEM (Number of Variables)
Lt=3*n        			    !MINIMUM LENGTH OF A MARKOW CHAIN
LtdOld=1.5*Lt  				!THIS VALUE IS USED ONLY IN THE FIRST DECREMENT OF THE CONTROL PARAMETER 
LtdOldP=LtdOld
NN=9*n     					!NUMBER OF CONFIGURATIONS

*DIM,Acost,ARRAY,NN   		!COST VALUE ARRAY, NN IS THE DIMENSION 
*DIM,order,ARRAY,NN   		!ORDER ARRAY, NN IS THE DIMENSION 
*DIM,Ah_rib,ARRAY,NN  	 	!INITIAL RIB HEIGHT ARRAY, NN IS THE DIMENSION 
*DIM,Aw_rib,ARRAY,NN  	 	!INITIAL RIB WIDTH ARRAY, NN IS THE DIMENSION 
*DIM,AL_rib,ARRAY,NN   		!INITIAL RIB LENGTH ARRAY, NN IS THE DIMENSION 
*DIM,Atheta_rib,ARRAY,NN   	!INTIAL ANGLE OF RIB ARRAY, NN IS THE DIMENSION
*DIM,ATsiWu,ARRAY,NN


ckini=200      				!INITIAL TEMPERATURE MODIFY FOR EVERY DIFFERENT MODEL
ck=ckini
alfamin=0.8  				!MINIMUM VALUE OF THE DECREMENT OF THE CONTROL PARAMETER
alfamax=0.9999  			!MAXIMUM VALUE OF THE DECREMENT OF THE CONTROL PARAMETER
epsin=0.00000002 			!STOPPING CRITERION (TEMPERATURE IS LOW ENOUGH)
epsinOut=0.0000002 			!STOPPING CRITERION (THE DIFFERENCE BETWEEN BEST AND WORST POINT IS NEGLIGIBLE)
NoIt=0         				!TOTAL NUMBER OF ITERATIONS EXCEPT GEOMETRICALLY UNACCEPTABLE MOVES
itotal=0       				!TOTAL NUMBER OF ITERATIONS 
mc=0           				!MARKOW CHAIN NUMBER
ann=-30        				!PARAMETER TO HANDLE VERY LOW NUMBERS

MassFirst=0					!INITIAL MASS
ObjectiveInitial=0			!INITIAL OBJECTIVE FUNCTION VALUE

!Defining Constraint Box for the Limits of the Rib
y4_min=-430e-3
y4_max=-20e-3
z4_min=-350e-3
z4_max=-20e-3

Hrib_min=0.020					!MINIMUM INITIAL NUMBER OF RIBS
Hrib_max=0.025					!MAXIMUM INITIAL NUMBER OF RIBS
Wrib_min=80e-3	 				!MINIMUM THICKNESS (WIDTH)
Wrib_max=160e-3     			!MAXIMUM THICKNESS (WIDTH)
L_rib_min=150e-3				!MINIMUM LENGTH
L_rib_max=270e-3				!MAXIMUM LENGTH 
Theta_rib_min=30				!MINIMUM ANGLE
Theta_rib_max=60				!MAXIMUM ANGLE

offset=-11e-4

rshrib_ini=(Hrib_max-Hrib_min)/3     	
rswrib_ini=(Wrib_max-Wrib_min)/3		
rsL_rib_ini=(L_rib_max-L_rib_min)/3		
rsthetarib_ini=(Theta_rib_max-Theta_rib_min)/3

rshrib=rshrib_ini			
rswrib=rswrib_ini			
rslrib=rsL_rib_ini			
rsthetarib=rsthetarib_ini

*END

!*******************************************************************************************************
!											INITIAL SHAPE
!*******************************************************************************************************

*CREATE,initialKey

!*******************************************************************************************************
!										SHAPE PROPERTIES OF THE BEAM
!*******************************************************************************************************

tc=2e-3			!thickness of web
ts=12e-3		!thickness of flange
H=500e-3		!Height of cross section 
Length=2000e-3	!Length of Beam
Wf=200e-3  		!Width of flange
LS=Length/10

!*******************************************************************************************************
!										SHAPE PROPERTIES OF THE RIB
!*******************************************************************************************************

H_rib=(Hrib_min+Hrib_max)/2			
W_rib=(Wrib_min+Wrib_max)/2			
L_rib=(L_rib_min+L_rib_max)/2			
Theta_rib=(Theta_rib_min+Theta_rib_max)/2

*END


!*******************************************************************************************************
!						DEFINITION OF THE KEYPOINTS AND GENERATION GEOMETRY
!*******************************************************************************************************

*CREATE,Definitionkey
/prep7

!*******************************************************************************************************
!											WEB AREA
!*******************************************************************************************************

L=Length
L1=Wf
L2=LS

K,1,-L1/2,0,(-L/2-L2)
K,2,0,0,(-L/2-L2)
K,3,L1/2,0,(-L/2-L2)
K,4,0,(-H+2*ts),(-L/2-L2)
K,5,-L1/2,(-H+2*ts),(-L/2-L2)
K,6,L1/2,(-H+2*ts),(-L/2-L2)

K,7,-L1/2,0,-L/2
K,8,0,0,-L/2
K,9,L1/2,0,-L/2
K,10,0,(-H+2*ts),-L/2
K,11,-L1/2,(-H+2*ts),-L/2
K,12,L1/2,(-H+2*ts),-L/2

K,13,-L1/2,0,0
K,14,0,0,0
K,15,L1/2,0,0
K,16,0,(-H+2*ts),0
K,17,-L1/2,(-H+2*ts),0
K,18,L1/2,(-H+2*ts),0

K,19,-L1/2,0,L/2
K,20,0,0,L/2
K,21,L1/2,0,L/2
K,22,0,(-H+2*ts),L/2
K,23,-L1/2,(-H+2*ts),L/2
K,24,L1/2,(-H+2*ts),L/2

K,25,-L1/2,0,(L/2+L2)
K,26,0,0,(L/2+L2)
K,27,L1/2,0,(L/2+L2)
K,28,0,(-H+2*ts),(L/2+L2)
K,29,-L1/2,(-H+2*ts),(L/2+L2)
K,30,L1/2,(-H+2*ts),(L/2+L2)

K,100,0,(-(H-2*ts)/2),0

*END

!*******************************************************************************************************
!					GENERATION OF THE AREAS AND CREATION OF THE LOCAL COORDINATE SYSTEMS
!*******************************************************************************************************

*CREATE,Areas

!---------------------------------------------------------------------------------------------------------
!									LOCAL COORDINATE SYSTEMS
!---------------------------------------------------------------------------------------------------------


*get,kp_max,KP,0,NUM,MAX
k,kp_max+1,0,0,-L/2
k,kp_max+2,0,0,(-L/2+1)
k,kp_max+3,L1/2,0,-L/2
k,kp_max+4,0,-(H-2*ts),-L/2
k,kp_max+5,0,-(H-2*ts),(-L/2+1)
k,kp_max+6,L1/2,-(H-2*ts),-L/2
k,kp_max+7,0,-(H-2*ts)/2,-L/2
k,kp_max+8,0,-(H-2*ts)/2,(-L/2+1)
k,kp_max+9,0,0,-L/2
k,kp_max+10,0,-(H-2*ts),L/2
k,kp_max+11,L1/2,-(H-2*ts),L/2
k,kp_max+12,0,-(H-2*ts)/2,L/2


CSKP,11,0,kp_max+1,kp_max+2,kp_max+3,    	!top flange coordinate
CSKP,12,0,kp_max+4,kp_max+5,kp_max+6, 		! bottom flange coordinate
CSKP,13,0,kp_max+7,kp_max+8,kp_max+9,    	! web coordinate
CSKP,14,0,kp_max+4,kp_max+6,kp_max+7,    	!edge at x=-l/2 
CSKP,15,0,kp_max+10,kp_max+11,kp_max+12,    !edge at x=L/2

!*******************************************************************************************************
!												WEB AREA
!*******************************************************************************************************

!Define Areas by Using Keyponts
csys,0

A,2,8,10,4
A,8,14,100,16,10
A,14,20,22,16,100
A,20,26,28,22

A,1,2,8,7
A,7,8,14,13
A,13,14,20,19
A,19,20,26,25
A,2,3,9,8
A,8,9,15,14
A,14,15,21,20
A,20,21,27,26

A,5,4,10,11
A,11,10,16,17
A,17,16,22,23
A,23,22,28,29
A,4,6,12,10
A,10,12,18,16
A,16,18,24,22
A,22,24,30,28

Asum
*get,Area0,AREA,0,AREA
!*******************************************************************************************************
!												Generate Ribs
!*******************************************************************************************************
csys,0
WPCSYS,-1,0 


rib1=0.20
rib2=0.40
rib3=0.20

k,151,offset,-(H-2*ts)/2+L_rib/2*sin(theta_rib)+W_rib/2*cos(theta_rib),L_rib/2*cos(theta_rib)-W_rib/2*sin(theta_rib)-rib1
k,152,H_rib,-(H-2*ts)/2+L_rib/2*sin(theta_rib),L_rib/2*cos(theta_rib)-rib1
k,153,offset,-(H-2*ts)/2+L_rib/2*sin(theta_rib)-W_rib/2*cos(theta_rib),L_rib/2*cos(theta_rib)+W_rib/2*sin(theta_rib)-rib1 

spline,151,152,153,,,,,1,-1,  ,-1,1

k,154,H_rib,-(H-2*ts)/2-L_rib/2*sin(theta_rib),-L_rib/2*cos(theta_rib)-rib1
l,152,154
adrag,51,52,,,,,53
    myerror1=_status
    *if,myerror1,ge,0.99,then
      reborder='yes'
    *endif
Ldele,53

k,159,offset,-(H-2*ts)/2+1.3*(L_rib/2+W_rib/2)*sin(theta_rib),1.3*(L_rib/2+W_rib/2)*cos(theta_rib)-rib1 
cskp,101,0,33,153,151,1,1
csys,101
kdele,154

spline,152,159,,,,,-1,0,0
adrag,51,52,,,,,53
    myerror2=_status
    *if,myerror2,ge,0.99,then
      reborder='yes'
    *endif

csys,0
k,162,offset,-(H-2*ts)/2-1.3*(L_rib/2+W_rib/2)*sin(theta_rib),-1.3*(L_rib/2+W_rib/2)*cos(theta_rib)-rib1
k,164,0,-(H-2*ts)/2-1.3*(L_rib/2+W_rib/2)*sin(theta_rib),-1.3*(L_rib/2+W_rib/2)*cos(theta_rib)-rib1
k,163,0,-(H-2*ts)/2+1.3*L_rib/2*sin(theta_rib),1.3*L_rib/2*cos(theta_rib)-rib1

csys,101
spline,32,162,,,,,1,0,0
adrag,54,57,,,,,64
    myerror3=_status
    *if,myerror3,ge,0.99,then
      reborder='yes'
    *endif

Ldele,64

KWPLAN,-1,     164,     163,     8
    myerror4=_status
    *if,myerror4,ge,0.99,then
      reborder='yes'
    *endif
	
block,-1.2*W_rib,1.6*L_rib+1.2*W_rib,-1.2*W_rib,1.2*W_rib,0,1.2*W_rib

asel,s,area,,21,26
asbv,all,1,,delete,delete
    myerror5=_status
    *if,myerror5,ge,0.99,then
      reborder='yes'
    *endif

allsel,all
aplot

al,88,89,92,85,83,91
    myerror6=_status
    *if,myerror6,ge,0.99,then
      reborder='yes'
    *endif
asba,2,21,,delete,delete
    myerror7=_status
    *if,myerror7,ge,0.99,then
      reborder='yes'
    *endif
allsel,all


csys,0
WPCSYS,-1,0 



agen,5,33,38,,,,-rib1,,0
    myerror8=_status
    *if,myerror8,ge,0.99,then
      reborder='yes'
    *endif


FLST,3,12,5,ORDE,6  
FITEM,3,2   
FITEM,3,21  
FITEM,3,23  
FITEM,3,-26 
FITEM,3,39  
FITEM,3,-44 
ARSYM,X,P51X, , , ,0,1  

al,60,62,65,57,54,64
asba,22,51,,delete,delete
    myerror9=_status
    *if,myerror9,ge,0.99,then
      reborder='yes'
    *endif
	
al,73,75,78,70,68,77
asba,52,22,,delete,delete
    myerror10=_status
    *if,myerror10,ge,0.99,then
      reborder='yes'
    *endif
	
al,98,100,103,95,93,102
asba,51,22,,delete,delete
    myerror11=_status
    *if,myerror11,ge,0.99,then
      reborder='yes'

    *endif

al,111,113,116,108,106,115
asba,52,22,,delete,keep
asba,1,22,,delete,delete
    myerror12=_status
    *if,myerror12,ge,0.99,then
      reborder='yes'
    *endif
allsel,all

FLST,3,30,5,ORDE,4  
FITEM,3,2   
FITEM,3,21  
FITEM,3,23  
FITEM,3,-50 
ARSYM,Z,P51X, , , ,0,0


    myerror13=_status
    *if,myerror13,ge,0.99,then
      reborder='yes'
    *endif


al,161,152,154,162,159,157
asba,3,81,,delete,delete

al,131,133,136,128,126,135
asba,82,3,,delete,delete

    myerror14=_status
    *if,myerror14,ge,0.99,then
      reborder='yes'
    *endif

al,144,146,148,149,141,139
asba,81,3,,delete,delete

al,170,172,175,167,165,174
asba,82,3,,delete,delete

    myerror15=_status
    *if,myerror15,ge,0.99,then
      reborder='yes'
    *endif
	
al,187,178,180,188,185,183
asba,81,3,,delete,keep
asba,4,3,,delete,delete
    myerror16=_status
    *if,myerror16,ge,0.99,then
      reborder='yes'
    *endif
allsel,all

aglue,all

    myerror17=_status
    *if,myerror17,ge,0.99,then
      reborder='yes'
    *endif
/Output,temp1,out,,Append
*VWRITE,'17  ',myerror17
(A4,F5.0)              
*VWRITE,jk,H_rib,W_rib,L_rib,Theta_rib,
(F5.0,F14.9,F14.9,F14.9,F14.9)  
 /Output         
*Cfclos

*END


!==============================================================!=========================================
!			GENERATION OF THE MESHES, APPLYING THE LOAD AND STABILIZATION BOUNDARY CONDITION
!==============================================================!=========================================
*CREATE,Meshes
/prep7
!---------------------------------------------------------------------------------------------------------
!		MESH
!---------------------------------------------------------------------------------------------------------

*do,i,5,12,1
CM,_Y,AREA  
ASEL, , , ,       i
CM,_Y1,AREA 
CMSEL,S,_Y  
!*  
CMSEL,S,_Y1 
AATT,       1, ,   1,       11,   1  
CMSEL,S,_Y  
CMDELE,_Y   
CMDELE,_Y1  
*enddo

!--------------------------------------------------------------------------------------------------------- 

*do,i,13,20,1
CM,_Y,AREA  
ASEL, , , ,       i
CM,_Y1,AREA 
CMSEL,S,_Y  
!*  
CMSEL,S,_Y1 
AATT,       1, ,   1,       12,   3  
CMSEL,S,_Y  
CMDELE,_Y   
CMDELE,_Y1  
*enddo

!--------------------------------------------------------------------------------------------------------- 

*do,i,1,4,1
CM,_Y,AREA  
ASEL, , , ,       i
CM,_Y1,AREA 
CMSEL,S,_Y  
!*  
CMSEL,S,_Y1 
AATT,       1, ,   1,       13,   2  
CMSEL,S,_Y  
CMDELE,_Y   
CMDELE,_Y1  
*enddo  

!---------------------------------------------------------------------------------------------------------
numcmp,area
*get,AreaTotalNo,area,0,count
*do,i,21,AreaTotalNo,1
CM,_Y,AREA  
ASEL, , , ,       i
CM,_Y1,AREA 
CMSEL,S,_Y  
!*  
CMSEL,S,_Y1 
AATT,       1, ,   1,       13,   2  
CMSEL,S,_Y  
CMDELE,_Y   
CMDELE,_Y1  
*enddo


!---------------------------------------------------------------------------------------------------------
!MESHING
!---------------------------------------------------------------------------------------------------------

Mesh_size=8e-3
Esize,8e-3				!Specifies the Default Number of Line Divisions
Amesh,All				!Generate Nodes and Area Elements within All Area by Using Defines Size in Esize Command
    myerror18=_status
    *if,myerror18,ge,2.99,then
      reborder='yes'
    *endif
/Output,temp1,out,,Append
*VWRITE,'18  ',myerror18
(A4,F5.0)              
/Output
*Cfclos

Asum
*get,TotalArea,AREA,0,AREA

*END


!*******************************************************************************************************
!											SOLUTION
!*******************************************************************************************************
*CREATE,Solution
/prep7
!*******************************************************************************************************
!								APPLIED LOAD AND BOUNDARY CONDITIONS
!*******************************************************************************************************

csys,0
!Defines DOF Constraints on Lines

DL,34,,UX,0
DL,34,,UY,0
DL,44,,UX,0
DL,44,,UY,0
DL,38,,UX,0
DL,38,,UY,0
DL,38,,UZ,0
DL,48,,UX,0
DL,48,,UY,0
DL,48,,UZ,0

DK,7,UX,0
DK,9,UX,0
DK,21,UX,0
DK,19,UX,0
DK,7,ROTZ,0
DK,9,ROTZ,0
DK,21,ROTZ,0
DK,19,ROTZ,0
!

!---------------------------------------------------------------------------------------------------------
!       				Master node for CERIG is created
!---------------------------------------------------------------------------------------------------------
K,5e6,0,0,0
ET,3,MASS21
R,1,0,0,0,0,0,0,
NUMSTR,NODE,5e6
KMESH,       5e6

!---------------------------------------------------------------------------------------------------------
!       				Rigid Region is created
!---------------------------------------------------------------------------------------------------------
NSEL,S,LOC,Z,(-(Mesh_size)/2),((Mesh_size)/2)
NSEL,R,LOC,Y,0,ts
NSEL,U, , ,    5e6 
CM,load_nodes,NODE,
CMSEL,S,load_nodes
NSEL,A, , ,    5e6 
CERIG,5e6,ALL,UY,RXYZ

ALLSEL

Force=-1   

F,5e6,FY,Force
Alls

 

FINISH  
/Solu
Alls
Antype,0							!Static Analysis
Solve
    myerror19=_status
    *if,myerror19,ge,1.99,then
      reborder='yes'
    *endif
/Output,temp1,out,,Append
*VWRITE,'19  ',myerror19
(A4,F5.0)              
/Output
*Cfclos

Finish

/solu
ANTYPE,1							!Linear Buckling Analysis
!*  
BUCOPT,LANB,6,0,0,CENTER
Pstres,on   
/STATUS,SOLU
SOLVE   
    myerror20=_status
    *if,myerror20,ge,1.99,then
      reborder='yes'
    *endif
/Output,temp1,out,,Append
*VWRITE,'20  ',myerror20
(A4,F5.0)              
/Output
*Cfclos

FINISH  

SAVE,'buckling','db'
FINISH

*END


!*******************************************************************************************************
!								CREATING THE MASS AND OBJECTIVE FUNCTIONS
!*******************************************************************************************************

*CREATE,costweightobj,mac


FINISH

/post1

!--------------------------------------------------------------------------------------------------------- 
!Mass0  =  Mass of I-beam without RIB
!Area0  =  Area of I-beam without RIB
!TotalArea = Total Area of I-beam With RIB

MassFirst=MassFirst+1
!Flange_AREA=L1*(Length+LS)*2
!Web_AREA=(H-2*ts)*(Length+LS)
!Rib_AREA=H_rib*((0.25*3.14*W_rib*W_rib)+(W_rib*L2))

Mass0=Density*9.81*Area0

!--------------------------------------------------------------------------------------------------------- 
!BuckleLF0 =  Buckling Load Factor of I-beam Without Rib
!BuckleLF  =  Buckling Load Factor of I-beam With Rib

BuckleLF0=5721.91

*GET,BuckleLF,MODE,1,FREQ

Mass=TotalArea*9.81*Density

*if,BuckleLF,le,0.0,then
   reborder='yes'
/Output,temp1,out,,Append
*VWRITE                                         
(T2,'Negative buckling factor',)                    
/Output
*Cfclos
*endif

ObjectiveInitial=ObjectiveInitial+1
Obj_F=0.85*(BuckleLF0/BuckleLF) + 0.15*(Mass0/Mass)

!---------------------------------------------------------------------------------------------------------

!FCTYP,ADD,TWSI
!nsort,fail,TWSI,0,0
!*get,MaxFail,sort,0,max

!---------------------------------------------------------------------------------------------------------
/OUTPUT,Objective_out,out,,APPEND
*VWRITE,jk,Obj_F,Mass0,Mass,BuckleLF,H_rib,W_rib,L_rib,Theta_rib,MaxFailTsaiWu
(F5.0,F14.6,F14.6,F14.6,F15.4,F12.7,F12.7,F12.7,F12.7,F9.5)
/OUTPUT
*CFCLOS

	/GRAPHICS,POWER
FINISH


*END

!*******************************************************************************************************
!					CALCULATING SHAPE PROPERTIES AND COST VALUES FOR THE CURRENT SHAPE
!*******************************************************************************************************

*CREATE,RefreshC

H_rib=Ah_rib(CurrentShape)			!INITIAL WAVE LENGTH
W_rib=Aw_rib(CurrentShape)			!INITIAL AMPLITUDE
L_rib=AL_rib(CurrentShape)			!FINAL AMPLITUDE
Theta_rib=Atheta_rib(CurrentShape)

*END

!*******************************************************************************************************
!							GENERATION OF RANDOM SHAPE PROPERTIES
!*******************************************************************************************************

*CREATE,RanPoint,mac

*DO,wwOut,1,10000

	*DO,ww,1,10000
		Step=rshrib*RAND(-1,1)
		TempH_rib=H_rib+Step
		*IF,TempH_rib,GT,Hrib_max,CYCLE
		*IF,TempH_rib,LT,Hrib_min,CYCLE
		*EXIT
	*ENDDO
	H_rib=TempH_rib

	!*IF,ww,GT,9998,THEN
	!	reborder='yes'
	!*ENDIF

	*DO,ww,1,10000
		Step=rswrib*RAND(-1,1)
		TempW_rib=W_rib+Step
		*IF,TempW_rib,GT,Wrib_max,CYCLE
		*IF,TempW_rib,LT,Wrib_min,CYCLE
		*EXIT
	*ENDDO
	W_rib=TempW_rib

	!*IF,ww,GT,9998,THEN
	!	reborder='yes'
	!*ENDIF
	
	*DO,ww,1,10000
		Step=rslrib*RAND(-1,1)
		TempL_rib=L_rib+Step
		*IF,TempL_rib,GT,L_rib_max,CYCLE
		*IF,TempL_rib,LT,L_rib_min,CYCLE
		*EXIT
	*ENDDO
	L_rib=TempL_rib

	!*IF,ww,GT,9998,THEN
	!	reborder='yes'
	!*ENDIF

	*DO,ww,1,10000
		Step=rsthetarib*RAND(-1,1)
		TempTheta_rib=Theta_rib+Step
		*IF,TempTheta_rib,GT,Theta_rib_max,CYCLE
		*IF,TempTheta_rib,LT,Theta_rib_min,CYCLE
		*EXIT
	*ENDDO
	Theta_rib=TempTheta_rib

	!*IF,ww,GT,9998,THEN
	!	reborder='yes'
	!*ENDIF
	
	
	Y4_rib=-(H-2*ts)/2+1.3*(L_rib/2+W_rib/2)*sin(Theta_rib)
	Z4_rib=1.3*(L_rib/2+W_rib/2)*cos(Theta_rib)-0.2


		*IF,Y4_rib,GT,y4_max,CYCLE
		*IF,Y4_rib,LT,y4_min,CYCLE
		*IF,Z4_rib,GT,z4_max,CYCLE
		*IF,Z4_rib,LT,z4_min,CYCLE
		*EXIT
	*ENDDO
	

	*IF,ww,GT,9998,THEN
		reborder='yes'
	*ENDIF

*END


!*******************************************************************************************************
!								FINDING THE BEST AND WORST POINTS
!*******************************************************************************************************

*CREATE,Reorder

order(1)=1
*DO,jj,2,NN
	order(jj)=jj
	*DO,ii,1,jj-1
		*IF,Acost(jj),LT,Acost(order(ii)),THEN
			*DO,kk,jj,ii+1,-1
				order(kk)=order(kk-1)
			*ENDDO
			order(ii)=jj
			*EXIT
		*ENDIF
	*ENDDO
*ENDDO
Best=order(1)
Worst=order(NN)
Worse=order(NN-3*n)
WorstBest=order(NN-3*n+1)

*END

!---------------------------------------------------------------------------------------------------------

/CLEAR,START

*DO,rknpp,1,1207
	Temp=RAND(0,1) 					!These lines prevent generating the same random moves 
*ENDDO 								!(ANSYS seems to Generate the same random points every time the program starts)

*USE,StartingPar 			!Starting Parameters
*USE,Initial 				!Initialization of the program
*USE,initialKey
*USE,Definitionkey

Ah_rib(1)=H_rib
Aw_rib(1)=W_rib
AL_rib(1)=L_rib
Atheta_rib(1)=Theta_rib


*USE,Areas 						!Generate area
*USE,Meshes 					!Generate the mesh
*USE,Solution	 				!Apply the boundary conditions and solve the problem

costweightobj

Acost(1)=Obj_F


/Output,tempo1,out,,Append
*VWRITE,jk,Obj_F,Mass,ww,jj
(F5.0,F14.6,F14.6,F5.0,F5.0)      
/Output
*Cfclos

jk=1
:reloop_j         !generating the other initial points
 jk=jk+1
 H_rib_old=H_rib
 W_rib_old=W_rib
 L_rib_old=L_rib
 Theta_rib_old=Theta_rib
 jj=1
 parsav,all
 *DOWHILE,jj
   /clear,start  
   /input,,parm  
   *use,Initial   				!initialization of the program 
   /prep7                     
   reborder='no'  
 	
   RanPoint
   RanPoint
	
   *USE,Definitionkey
   *IF,reborder,EQ,'no',THEN 
     *USE,Areas
   *ENDIF
   
    
   *IF,reborder,EQ,'no',THEN 
 	 *USE,Meshes
   *ENDIF
   *IF,reborder,EQ,'no',THEN
	 *USE,Solution
   *ENDIF
   *IF,reborder,EQ,'no',THEN
	 costweightobj
   *ENDIF
   *IF,reborder,EQ,'no',THEN
 	 jj=0
   *ELSE
	 H_rib=H_rib_old
	 W_rib=W_rib_old
	 L_rib=L_rib_old
	 Theta_rib=Theta_rib_old
   *ENDIF
 *enddo
   	   
 Ah_rib(jk)=H_rib
 Aw_rib(jk)=W_rib
 AL_rib(jk)=L_rib
 Atheta_rib(jk)=Theta_rib
 Acost(jk)=Obj_F
 ATsiWu(jk)=MaxFailTsaiWu
	
/Output,tempo2,out,,Append
*VWRITE,jk,Obj_F,H_rib,W_rib,L_rib,Theta_rib,ww,jj,Mass,MaxFailTsaiWu
(F5.0,F14.6,F14.6,F14.6,F14.6,F14.6,F5.0,F5.0,F12.6,F12.6)              
/Output
*Cfclos
  
*if,jk,lt,NN,:reloop_j  

*USE,Reorder   !Finding the best and worst point

exitLoop='true' 
exitOut='false'
mc=0
PARSAV,ALL
FINISH

/Output,BestConfigurations,out,,Append
*VWRITE
(T2,'mc    BestObj BestBuckleLF BestMass  BestTW   Best_H   Best_w   Best_L  Best_theta')            
/Output
*Cfclos
 
/COPY,,parm,,inipoint,,,


!---------------------------------------------------------------------------------------------------------
